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ABSTRACT 


This  note  provides  a  communication-theoretic  rationale  for  a 
successful  nonlinear  processing  scheme  due  to  Shimshoni  and  Smith, 
by  showing  its  resemblance  to  minimum -rms -err or  estimation.  The 
latter  filtering,  although  "optimum"  under  certain  assumptions,  is 
rather  more  complicated  in  its  implementation. 

As  in  all  such  filtering  studies,  the  mean-square  error  is  made 
up  of  bias,  or  signal  distortion,  and  variance  (which  in  linear  filtering 
is  attributable  only  to  the  noise).  It  would  be  interesting  to  see  whether 
the  filtering  developed  here  could  be  modified  to  allow  a  controlled 
tradeoff  between  these  two  sources  of  error. 
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Shimshoni  and  Smith  have  described  a  procedure  that  takes  advantage  of  the 
linear  polarization  in  a  P-wave  arrival  to  enhance  its  visibility  against  noise  that  is 
assumed  to  be  unpolarized.  Their  technique  is  to  multiply  the  output  of  the  vertical 
seismometer  by  a  running  history  of  the  short-term,  zero-lag  crosscorrelation  taken 
between  the  vertical  and  horizontal  seismometer  outputs,  with  the  horizontal  seis¬ 
mometer  aligned  to  the  radial  component  of  P-wave  particle  motion. 

The  purpose  of  his  note  is  to  demonstrate  that  while  the  Shimshoni -Smith 
procedure  may  at  first  seem  quite  ad  hoc,  it  is  in  fact  ’’cousin”  to  a  nonlinear  pro¬ 
cessing  scheme  which,  under  assumptions  to  be  given,  is  optimal  in  the  sense  that  its 
output  is  (proportional  to)  the  minimum -variance  estimate  of  the  P-wave  signature. 
Although  it  may  not  be  worthwhile  to  implement  this  optimum  scheme,  its  similarity 
to  the  Shimshoni -Smith  procedure  at  least  can  serve  to  explain  the  latter's  marked 
success  and  support  its  philosophy. 

We  model  the  problem  as  follows.  The  outputs  of  the  vertical  and  horizontal 
seismometers  are  assumed  to  contain  independent  gaussian  noises  (i.e.  ,  the  noise  in 
the  particle  motion  is  unpolarized)  that  are  white  and  of  equal  intensity  (these  restric¬ 
tions  are  not  essential,  but  are  made  merely  for  analytical  convenience  in  this  cursory 
study).  The  P-wave  is  taken  to  be  the  result  of  translating  a  sample  of  non -stationary 
gaussian  noise,  which  initially  vanishes  outside  an  interval  (0,A)  and  has  a  known  coi 
relation  function  9p(t,T)  inside,  by  an  unknown  amount  tg,  so  that  the  (unknown -waveform) 
P-wave  signature  occupies  the  interval  (ts,tg+A).  Its  starting  time  tg  is  assumed  to 
have  a  known  a  priori  probability  density  over  some  finite  interval  (t0,t0+T),  which 
vanishes  outside.  The  vertical  and  horizontal  seismometers  are  constructed  to 
respond  identically  to  the  arriving  P-wave  except  for  sine  and  cosine  gain  factors  that 
depend  on  the  arrival  angle.  To  form  a  minimum -variance  estimate  of  the  P-wave 
signature,  one  is  forced,  as  in  Wiener  linear  filtering,  to  assume  some  a  priori 
statistics  for  unknown  quantities  as  above;  if  the  arrival  angle  is  unknown  as  well,  we 
must  likewise  assign  it  an  a  priori  distribution. 
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At  least  in  the  important  case  of  discerning  the  pP  arrival,  it  appears  safe  to 
assume  that  the  arrival  angle  is  known,  this  information  being  provided  by  the  generally 
stronger  main-P  arrival.  We  shall  depart  from  the  Shimshoni-Smith  context  by 
assuming  the  arrival  angle  to  be  known.  If  it  is  indeed  unknown,  one  would  obtain  the 
minimum -variance  estimate  by  applying  the  following  analysis  to  each  of  a  continuum 
of  hypothesized  arrival  angles,  and  then  averaging  the  resulting  infinitely-many 
minimum -variance  estimates  (each  conditioned  on  a  different  hypothesized  angle)  over 
the  a  priori  distribution  of  angles. 

At  the  outset,  we  know  that  the  minimum -variance  pP-wave  signature  estimator 
must  be  nonlinear,  for  the  signature  has  non-gaussian  statistics  as  a  result  of  the 
averaging  of  its  conditionally  gaussian  (on  ts)  statistics  over  the  random  starting  time 
ts.  Hence  there  should  be  no  prejudice  a  priori  against  procedures  like  those  of 
Shimshoni-Smith.  However,  we  want  to  see  where  we  are  led  through  a  synthesis  of 
the  nonlinear  processing  from  basic  statistical  principles  (admitting  that  this  requires 
questionable  assumptions  of  a  priori  statistics),  avoiding  ad  hoc  notions  even  though 
the  outcome  of  this  abstract  synthesis  may  have  a  straightforward  heuristic  interpre¬ 
tation. 

To  begin,  we  are  given  the  respective  outputs  Wy(t)  and  Wjj(t)  of  the  vertical 
and  horizontal  seismometers,  observed  over  all  time  (this  assumption  is  not  essential). 
Knowing  the  pP  arrival  angle,  we  may  form  a  new  pair  of  waveforms  w(t)  and  w(t) 
through  a  (time -independent)  rotational  transformation  of  Wy(t)  and  Wjj(t),  in  such  a 
way  that  w(t)  contains  only  noise  (independent  of  that  in  w(t)  by  virtue  of  the  assumed 
equal  noise  intensities  for  wy(t)  and  wpj(t)  and  the  fact  that  the  transformation  is  a 
rotation)  and  no  pP-wave  energy.  (In  effect,  we  construct  a  new  pair  of  orthogonally 
aligned  seismometers,  the  one  whose  output  is  w(t)  responding  only  at  right-angles  to 
the  direction  of  signal  particle  motion.)  Nothing  is  lost  in  performing  such  a  rotation, 
for  it  is  a  reversible  operation;  therefore,  the  minimum  P-signature  error  variance 
obtainable  by  observing  w(t)  and  w{t)  is  identical  to  that  obtainable  by  observing  the 
original  outputs  wy(t)  and  wj^(t). 
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It  is  a  universal  rule,  irrespective  of  whether  or  not  the  statistics  that  are 
involved  are  gaussian,  that  the  minimum -variance  estimator  (which  generally  must 
suffer  some  bias  or  signal  distortion,  as  in  Wiener  estimation)  always  computes  the 
conditional  mean  of  the  signal  whose  estimate  is  sought,  where  the  conditioning  is  on 
whatever  observations  are  available.  Thus  our  optimum  estimator  is 

CD 

xW,W(t)  =  x(t)'w  W  =  J’  xtP[xt/W,W]  (1) 

where  x(t)  is  the  pP-wave  signature  whose  minimum -variance  estimate  is  sought  from 
the  observations  of  w(t)  and  w(t),  and  this  latter  pair  of  histories  is  for  conciseness 
denoted  by  W  and  W,  respectively.  (In  averaging,  the  notation  Xj.  rather  than  x(t)  is 
used  to  indicate  that  in  this  instance  we  are  concerned  with  only  a  single  instant  of 
time. ) 

We  now  have  recourse  to  the  result 


p[xt/WfW] 


EL! y.»^»xt)  _  Ptw]  p[w,xt]  _ 

p[W,W]  plW]  p[W]  ^ 


(2) 


Philosophically,  the  fact  that  w(t)  contains  only  noise  and  that  this  noise  is  independent 
of  that  in  w(t)  causes  it  to  play  no  role  in  the  estimation  of  the  P-wuve  signature,  as 
has  just  been  demonstrated  mathematically. 

Using  the  identity: 


P[xt/W] 


t+T 


p[x./W,ts]  p[tg/W]  dtg 


(3) 


and  Bayes's  rule: 


P(ts/W) 


p(W/tg)  p(tg) 

p(W) 


(4) 
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our  estimate  becomes 


_  -  tQ+T  00 

*w  w(t)  =  *W(t*  =  lp(W>]  J*  {  J*  *  P[xt/W,tsj  dxj  p(W/ts)  p(tg)  dtg  (5) 

’  fo  — 

The  braced  inner  integral  in  (5),  which  may  conveniently  be  denoted  by  x^  (t),  is 
identical  to  the  linear  Wiener  estimate  of  x(t),  conditioned  on  both  the  observation  of 
w(t)  and  the  hypothesis  that  the  starting  time  for  x(t)  is  tg,  for  it  can  be  shown  that 
plXj/Wjts]  is  gaussian  in  x  and  thence  that 

00 

xw  t  (t)  =  J  w(t)  h(T-ts,t-ts)  dT  (6) 

9  S  —  oo 

where  the  (symmetric,  nonrealizable)  linear  operator  h(T,t)  solves  the  integral 
equation 

A 

J*  (cpp(t,r)  +  (No/2)  6(t-x)]  h(T,a)dx  =  cpp(T,a)  ;  O^t.a^A  (7) 

0 

and  vanishes  for  either  Tort  outside  (0,A).  Here  cpp(t,  t)  is  the  correlation  function  of 
the  gaussian  P-wave  signature  for  zero  translation  (tg  =  0),  5(t)  is  the  Dirac  delta - 
function  or  unit  impulse,  and  N0  is  the  (single -sided,  physical)  spectral  density  of  the 
white  noise  in  the  observation  w(t). 

A 

It  is  now  clear  from  (5)  that,  taking  into  account  the  dependence  of  on 

A 

p(W/tJ  as  well  as  (linearly)  on  w(t)  through  x,„  (t),  the  minimum -variance  estimator 

b  w»ls  -1 

will  be  nonlinear  in  w(t).  (Henceforth,  we  shall  ignore  the  outside  factor  [p(W)  ]  in 

(5),  as  for  any  given  observation  it  represents  a  fixed  gain.  Thus  our  estimate  will 

now  be  proportional  to  the  minimum -variance  estimate,  rather  than  strictly  equaling 

it- -additional  computation  could  supply  p(W)  if  it  were  needed,  but  ordinary  AGC  along 

with  compensation  by  the  human  eye  should  suffice  unless  the  magnitude  of  pP  must  be 

estimated. ) 
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It  remains  to  determine  p(W/tg),  the  a  priori  probability  of  obtaining  the  given 
observation  under  the  hypothesis  that  the  (conditionally)  gaussian  P-wave  signature 
commences  at  tg.  Applying  the  results  for  the  detection  of  a  gaussian  signal  in  white 
gaussian  noise  (e.g.  ,  Price,  IRE  Trans,  on  Info.  Theory,  December  1956),  we  have 

00  00 

p(W/tg)  ~  exp  {  J  w(t)  [  j  w(T)  h(T-tg,t-ts)  dT  ]  dt  }  (8) 

O  —  CD  —  00 

where  again  we  are  not  concerned  with  the  (here  w(t)-independent)  constant  of  propor¬ 
tionality.  Assembling  (5)-(8),  we  arrive  at  the  block  diagram  (page  7)  for  the  pro¬ 
cessor  whose  output  is  proportional  to  the  minimum -variance  estimate  of  the  P-wave, 
and  which  exhibits  linear  and  quadratic  operations,  and  a  "bootstrap"  connection 
between  them  not  unlike  that  of  Shimshoni -Smith. 

With  regard  to  the  quadratic  operation,  Shimshoni -Smith  multiply  the  original 
seismometer  outputs  Wy(t)  and  wj.j(t)  and  integrate  the  product  for  a  short  period. 
Recalling  that  the  input  to  the  processor  shown  in  the  figure  results  from  a  rotation 
and  is  thus  of  the  form  w(t)  =  awy(t)  +  bwjj(t),  we  note  that  such  a  product  (with  some 
filtering)  is  among  the  components  appearing  at  the  output  of  the  uppermost  multiplier, 
and  that  for  any  given  t  the  indicated  infinite-time  integration  is  in  fact  only  finite  by 
virtue  of  the  A-limited  filter  memories. 
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